rm(list = ls())
suppressPackageStartupMessages(library(ggplot2))
suppressPackageStartupMessages(library(tidyverse))
suppressPackageStartupMessages(library(plyr))
suppressPackageStartupMessages(library(dplyr))
suppressPackageStartupMessages(library(knitr))
suppressPackageStartupMessages(library(kableExtra))
suppressPackageStartupMessages(library(ggpubr))
suppressPackageStartupMessages(library(pander))
suppressPackageStartupMessages(library(jtools))
# suppressPackageStartupMessages(library(stringr))
# suppressPackageStartupMessages(library(broom.mixed))
# suppressPackageStartupMessages(library(sjPlot))
# suppressPackageStartupMessages(library(huxtable))
# suppressPackageStartupMessages(library(ggstance))pd <- position_dodge(0.3) # move them .05 to the left and right
# Multiple plot function
#
# ggplot objects can be passed in ..., or to plotlist (as a list of ggplot objects)
# - cols: Number of columns in layout
# - layout: A matrix specifying the layout. If present, 'cols' is ignored.
#
# If the layout is something like matrix(c(1,2,3,3), nrow=2, byrow=TRUE),
# then plot 1 will go in the upper left, 2 will go in the upper right, and
# 3 will go all the way across the bottom.
#
multiplot <- function(..., plotlist=NULL, file, cols=1, layout=NULL) {
library(grid)
# Make a list from the ... arguments and plotlist
plots <- c(list(...), plotlist)
numPlots = length(plots)
# If layout is NULL, then use 'cols' to determine layout
if (is.null(layout)) {
# Make the panel
# ncol: Number of columns of plots
# nrow: Number of rows needed, calculated from # of cols
layout <- matrix(seq(1, cols * ceiling(numPlots/cols)),
ncol = cols, nrow = ceiling(numPlots/cols))
}
if (numPlots==1) {
print(plots[[1]])
} else {
# Set up the page
grid.newpage()
pushViewport(viewport(layout = grid.layout(nrow(layout), ncol(layout))))
# Make each plot, in the correct location
for (i in 1:numPlots) {
# Get the i,j matrix positions of the regions that contain this subplot
matchidx <- as.data.frame(which(layout == i, arr.ind = TRUE))
print(plots[[i]], vp = viewport(layout.pos.row = matchidx$row,
layout.pos.col = matchidx$col))
}
}
}
numericcharacters <- function(x) {
!any(is.na(suppressWarnings(as.numeric(x)))) & is.character(x)
}
scale1 <- function(x) scale(x)[,1]
construct_coef_df = function(dvs,df){
fhab_mds = lapply(dvs, function(x) {
lm(substitute(i ~ fhab + study_t + scanner_t, list(i = as.name(x))), data = df)
})
fhab_sums = lapply(fhab_mds, summary)
femo_mds = lapply(dvs, function(x) {
lm(substitute(i ~ femo + study_t + scanner_t, list(i = as.name(x))), data = df)
})
femo_sums = lapply(femo_mds, summary)
fhab_betas = c()
femo_betas = c()
fhab_se = c()
femo_se = c()
for (i in c(1:length(dvs))){
fhab_betas = c(fhab_betas,fhab_sums[[i]]$coefficients[2,1])
femo_betas = c(femo_betas,femo_sums[[i]]$coefficients[2,1] )
fhab_se = c(fhab_se, fhab_sums[[i]]$coefficients[2,2])
femo_se = c(femo_se, femo_sums[[i]]$coefficients[2,2])
}
coef_df = data.frame("IV" = c(rep("Habitual FB use", length(dvs)), rep("Emotional FB use", length(dvs))),
"DV" = c(dvs,dvs),
"beta" = c(fhab_betas, femo_betas),
"se" = c(fhab_se, femo_se))
coef_df = coef_df %>%
mutate(IV = factor(coef_df$IV, levels=c("Habitual FB use", "Emotional FB use")))
return(coef_df)
}
plot_coef_df = function(coef_dv){
g1 = ggplot(coef_df, aes(x=DV, y=beta, colour=IV)) +
geom_errorbar(aes(ymin= beta-1.96*se, ymax= beta+1.96*se), width=.1, position=pd) +
geom_point(position=pd) + theme_pubr() +
scale_color_manual(values = c("Habitual FB use" = "dodgerblue", "Emotional FB use" = "forestgreen")) +
geom_hline(yintercept=0, linetype="dashed",
color = "grey", size=1) +
rotate_x_text(45)
return(g1)
}
construct_coef_df_iri = function(dvs,df){
mds1 = lapply(dvs, function(x) {
lm(substitute(i ~ IRI_Perspective_Taking + study_t + scanner_t, list(i = as.name(x))), data = df)
})
sums1 = lapply(mds1, summary)
mds2 = lapply(dvs, function(x) {
lm(substitute(i ~ IRI_Fantasy + study_t + scanner_t, list(i = as.name(x))), data = df)
})
sums2 = lapply(mds2, summary)
mds3 = lapply(dvs, function(x) {
lm(substitute(i ~ IRI_Empathic_Concern + study_t + scanner_t, list(i = as.name(x))), data = df)
})
sums3 = lapply(mds3, summary)
mds4 = lapply(dvs, function(x) {
lm(substitute(i ~ IRI_Personal_Distress + study_t + scanner_t, list(i = as.name(x))), data = df)
})
sums4 = lapply(mds4, summary)
betas1 = c()
betas2 = c()
betas3 = c()
betas4 = c()
se1 = c()
se2 = c()
se3 = c()
se4 = c()
for (i in c(1:length(dvs))){
betas1 = c(betas1,sums1[[i]]$coefficients[2,1])
betas2 = c(betas2,sums2[[i]]$coefficients[2,1])
betas3 = c(betas3,sums3[[i]]$coefficients[2,1])
betas4 = c(betas4,sums4[[i]]$coefficients[2,1])
se1 = c(se1, sums1[[i]]$coefficients[2,2])
se2 = c(se2, sums2[[i]]$coefficients[2,2])
se3 = c(se3, sums3[[i]]$coefficients[2,2])
se4 = c(se4, sums4[[i]]$coefficients[2,2])
}
coef_df = data.frame("IV" = c(rep("Perspective taking", length(dvs)),
rep("Fantasy", length(dvs)),
rep("Empathic concern", length(dvs)),
rep("Personal distress", length(dvs))),
"DV" = c(dvs,dvs,dvs,dvs),
"beta" = c(betas1, betas2,betas3,betas4),
"se" = c(se1, se2,se3,se4))
coef_df = coef_df
return(coef_df)
}
plot_coef_iri = function(coef_dv){
g1 = ggplot(coef_df, aes(x=DV, y=beta, colour=IV)) +
geom_errorbar(aes(ymin= beta-1.96*se, ymax= beta+1.96*se), width=.1, position=pd) +
geom_point(position=pd) + theme_pubr() +
geom_hline(yintercept=0, linetype="dashed",
color = "grey", size=1) +
rotate_x_text(45)
return(g1)
}
corstar <- function(x, y = NULL, use = "pairwise", method = "pearson", round = 3, row.labels, col.labels, ...) {
require(psych)
ct <- corr.test(x, y, use, method) # calculate correlation
r <- ct$r # get correlation coefs
p <- ct$p # get p-values
stars <- ifelse(p < .001, "***", ifelse(p < .01, "** ", ifelse(p < .05, "* ", " "))) # generate significance stars
m <- matrix(NA, nrow = nrow(r) * 2, ncol = ncol(r) + 1) # create empty matrix
rlab <- if(missing(row.labels)) rownames(r) else row.labels # add row labels
clab <- if(missing(col.labels)) {
if(is.null(colnames(r)))
deparse(substitute(y))
else
colnames(r)
} else {
col.labels # add column labels
}
rows <- 1:nrow(m) # row indices
cols <- 2:ncol(m) # column indices
odd <- rows %% 2 == 1 # odd rows
even <- rows %% 2 == 0 # even rows
m[odd, 1] <- rlab # add variable names
m[even, 1] <- rep("", sum(even)) # add blank
m[odd, cols] <- paste(format(round(r, round), nsmall = round, ...), stars, sep = "") # add r coefs
m[even, cols] <- paste("(", format(round(p, round), nsmall = round, ...), ")", sep = "") # add p-values
colnames(m) <- c(" ", clab) # add colnames
m # return matrix
}path = "~/Box Sync/CurrentProjects_Penn/TPS12_BART/99_fomo/github/"
setwd(path)
## socialMedia data
socialMedia = read.csv(paste0(path, '0_facebook/tps2_facebook.csv'),na.strings=c("","NA"), stringsAsFactors = FALSE)[-1,]
socialMedia = socialMedia %>%
mutate(TPS_ID = as.factor(TPS_ID)) %>%
mutate_if(is.character,as.numeric) %>%
mutate(fhab = rowMeans(.[,10:19], na.rm = FALSE),
femo = rowMeans(.[c(8,9,20:39)], na.rm = FALSE), ## include the two "feel connected" questions in Emotional FB use.
TPS_ID = as.character(TPS_ID)) %>%
select(TPS_ID, fhab, femo)
cyber = read.csv(paste0(path, '0_cyberball/cyberball_roi.csv'), stringsAsFactors = FALSE)
age = read.csv(paste0(path, '0_surveys/tps12_scanner_study_age.csv'), stringsAsFactors = FALSE) %>%
select(TPS_ID, age, study_t, scanner_t)
iri = read.csv(paste0(path, '0_surveys/iri_tps12.csv'), stringsAsFactors = FALSE)
iri = iri %>%
mutate(IRI_mean = rowMeans(select(iri, starts_with("IRI_")), na.rm = FALSE))
nts = read.csv(paste0(path, '0_surveys/nts.csv'), stringsAsFactors = FALSE)
df = cyber %>%
left_join(socialMedia, by = "TPS_ID") %>%
left_join(age, by = "TPS_ID") %>%
left_join(iri, by = "TPS_ID") %>%
left_join(nts, by = "TPS_ID")
# Remove NAs
df = df %>% drop_na(age,fhab,femo, NTS_Belongingness,
NTS_SelfEsteem, NTS_Control,
NTS_MeaningfulExistence, study_t, scanner_t,
IRI_Perspective_Taking, IRI_Fantasy,
IRI_Empathic_Concern, IRI_Personal_Distress)We include 59 participants in this report.
1 - Strongly disagree 7 - Strongly agree
I get an urge to post on Facebook from my computer or phone the moment when … 3. something makes me feel amused
4. something makes me feel surprised 5. something makes me feel awed 6. something makes me feel loved 7. something makes me feel proud
8. something makes me feel excited
9. something makes me feel grateful 10. something makes me feel happy
11. something makes me feel inspired
12. something makes me feel confident
13. something makes me feel angry
14. something makes me feel nervous 15. something makes me feel awkward 16. something makes me feel stressed
17. something makes me feel jealous 18. something makes me feel lonely
19. something makes me feel scared
20. something makes me feel upset
21. something makes me feel ashamed 22. something makes me feel guilty
Scrit to obtain the ROI values at jupyter hub.
roi_names = colnames(df)[c(3:20,23)]
source = c("functional ROI",
"functional ROI",
"functional ROI",
"functional ROI",
"Dufour et al., 2013","Dufour et al., 2013","Dufour et al., 2013",
"Dufour et al., 2013","Dufour et al., 2013","Dufour et al., 2013",
"Dufour et al., 2013",
"Vijayakumar et al., 2017","Vijayakumar et al., 2017",
"Vijayakumar et al., 2017","Vijayakumar et al., 2017",
"Vijayakumar et al., 2017","Vijayakumar et al., 2017",
"Vijayakumar et al., 2017","Vijayakumar et al., 2017")
img_folder = '/Users/Rui/Box Sync/CurrentProjects_Penn/TPS12_BART/99_fomo/data/cyberball/ROIs/png/'
imgs = paste0(img_folder, roi_names, '.png')
rois = data.frame("name" = roi_names, "source" = source, "image" = imgs)
rois %>%
mutate(
image %>% pander::pandoc.image.return()
) %>%
select(-image) %>%
pander()| name | source |
|---|---|
| subacc_roi_mask | functional ROI |
| rInsula_roi_mask | functional ROI |
| lInsula_roi_mask | functional ROI |
| functional_combined_roi_mask | functional ROI |
| dmpfc_roi_mask | Dufour et al., 2013 |
| mmpfc_roi_mask | Dufour et al., 2013 |
| vmpfc_roi_mask | Dufour et al., 2013 |
| mpfc_roi_mask | Dufour et al., 2013 |
| precuneus_roi_mask | Dufour et al., 2013 |
| lTpj_roi_mask | Dufour et al., 2013 |
| rTpj_roi_mask | Dufour et al., 2013 |
| rSts_roi_mask | Vijayakumar et al., 2017 |
| saxe_combined_mask | Vijayakumar et al., 2017 |
| ns_ment_mask | Vijayakumar et al., 2017 |
| ns_ment_clust1_mask | Vijayakumar et al., 2017 |
| ns_ment_clust2_mask | Vijayakumar et al., 2017 |
| ns_ment_clust3_mask | Vijayakumar et al., 2017 |
| ns_ment_clust4_mask | Vijayakumar et al., 2017 |
| ns_ment_clust7_mask | Vijayakumar et al., 2017 |
| image %>% pander::pandoc.image.return() |
|---|
Perspective taking
ggdensity(df, x = "IRI_Perspective_Taking",
add = "mean", rug = TRUE,
color = "study_t", fill = "study_t")Perspective taking
ggdensity(df, x = "IRI_Personal_Distress",
add = "mean", rug = TRUE,
color = "study_t", fill = "study_t")Empathic concern
ggdensity(df, x = "IRI_Empathic_Concern",
add = "mean", rug = TRUE,
color = "study_t", fill = "study_t")Fantasy
Lower scores indicate greater distress
Belongingness
ggdensity(df, x = "NTS_Belongingness",
add = "mean", rug = TRUE,
color = "study_t", fill = "study_t")Self esteem
Control
MeaningfulExistence
ggdensity(df, x = "NTS_MeaningfulExistence",
add = "mean", rug = TRUE,
color = "study_t", fill = "study_t")Stars indicate signiciant correlation, and (x) indicate p values.
m = corstar(df[,c("IRI_Perspective_Taking","IRI_Fantasy",
"IRI_Empathic_Concern","IRI_Personal_Distress",
"NTS_Belongingness","NTS_SelfEsteem",
"NTS_Control","NTS_MeaningfulExistence")])
m[c(1,2),2] = ""
m[c(1:4),3] = ""
m[c(1:6),4] = ""
m[c(1:8),5] = ""
m[c(1:10),6] = ""
m[c(1:12),7] = ""
m[c(1:14),8] = ""
df_m = as.data.frame(m[,1:8])
rownames(df_m) = NULL
kable(df_m) %>%
kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive"))| IRI_Perspective_Taking | IRI_Fantasy | IRI_Empathic_Concern | IRI_Personal_Distress | NTS_Belongingness | NTS_SelfEsteem | NTS_Control | |
|---|---|---|---|---|---|---|---|
| IRI_Perspective_Taking | |||||||
| IRI_Fantasy | 0.220 | ||||||
| (0.094) | |||||||
| IRI_Empathic_Concern | 0.495*** | 0.311* | |||||
| (0.000) | (0.017) | ||||||
| IRI_Personal_Distress | -0.144 | 0.115 | 0.202 | ||||
| (0.275) | (0.386) | (0.124) | |||||
| NTS_Belongingness | 0.246 | 0.045 | 0.114 | -0.181 | |||
| (0.060) | (0.736) | (0.391) | (0.171) | ||||
| NTS_SelfEsteem | 0.068 | 0.202 | -0.036 | -0.217 | 0.486*** | ||
| (0.610) | (0.124) | (0.789) | (0.099) | (0.000) | |||
| NTS_Control | 0.086 | 0.166 | -0.041 | -0.058 | 0.680*** | 0.612*** | |
| (0.519) | (0.210) | (0.756) | (0.665) | (0.000) | (0.000) | ||
| NTS_MeaningfulExistence | 0.287* | 0.273* | 0.140 | -0.046 | 0.728*** | 0.438*** | 0.709*** |
| (0.028) | (0.037) | (0.291) | (0.727) | (0.000) | (0.001) | (0.000) |
dvs = c("subacc_roi_mask","rInsula_roi_mask",
"lInsula_roi_mask", "functional_combined_roi_mask")
coef_df = construct_coef_df(dvs, df) %>%
mutate(DV = gsub("_roi_mask", "", DV))
plot_coef_df(coef_dv)dvs = c("vmpfc_roi_mask","dmpfc_roi_mask","mmpfc_roi_mask", "precuneus_roi_mask", "lTpj_roi_mask","rTpj_roi_mask" ,"rSts_roi_mask",
"saxe_combined_mask")
coef_df = construct_coef_df(dvs, df) %>%
mutate(DV = gsub("_roi_mask", "", DV))
plot_coef_df(coef_dv)df = df %>%
dplyr::rename(.,neurosynth_mask = ns_ment_mask,
ns_RSTS = ns_ment_clust1_mask,
ns_RTPJ = ns_ment_clust2_mask,
ns_DMPFC = ns_ment_clust3_mask,
ns_PCC = ns_ment_clust4_mask,
ns_VMPFC = ns_ment_clust5_mask,
ns_hippocampus = ns_ment_clust6_mask,
ns_LTPJ = ns_ment_clust7_mask,
ns_LSTS = ns_ment_clust8_mask)
dvs = c("neurosynth_mask","ns_RSTS", "ns_RTPJ",
"ns_DMPFC","ns_PCC" ,"ns_VMPFC",
"ns_hippocampus","ns_LTPJ","ns_LSTS")
coef_df = construct_coef_df(dvs, df) %>%
mutate(DV = gsub("_mask", "", DV))
plot_coef_df(coef_dv)dvs = c("vmpfc_roi_mask","dmpfc_roi_mask","mmpfc_roi_mask", "precuneus_roi_mask", "lTpj_roi_mask","rTpj_roi_mask" ,"rSts_roi_mask",
"saxe_combined_mask")
coef_df = construct_coef_df_iri(dvs, df) %>%
mutate(DV = gsub("_roi_mask", "", DV))
plot_coef_iri(coef_dv)## function
construct_coef_df_nts = function(dvs,df){
mds1 = lapply(dvs, function(x) {
lm(substitute(i ~ NTS_Belongingness + study_t + scanner_t, list(i = as.name(x))), data = df)
})
sums1 = lapply(mds1, summary)
mds2 = lapply(dvs, function(x) {
lm(substitute(i ~ NTS_SelfEsteem + study_t + scanner_t, list(i = as.name(x))), data = df)
})
sums2 = lapply(mds2, summary)
mds3 = lapply(dvs, function(x) {
lm(substitute(i ~ NTS_Control + study_t + scanner_t, list(i = as.name(x))), data = df)
})
sums3 = lapply(mds3, summary)
mds4 = lapply(dvs, function(x) {
lm(substitute(i ~ NTS_MeaningfulExistence + study_t + scanner_t, list(i = as.name(x))), data = df)
})
sums4 = lapply(mds4, summary)
betas1 = c()
betas2 = c()
betas3 = c()
betas4 = c()
se1 = c()
se2 = c()
se3 = c()
se4 = c()
for (i in c(1:length(dvs))){
betas1 = c(betas1,sums1[[i]]$coefficients[2,1])
betas2 = c(betas2,sums2[[i]]$coefficients[2,1])
betas3 = c(betas3,sums3[[i]]$coefficients[2,1])
betas4 = c(betas4,sums4[[i]]$coefficients[2,1])
se1 = c(se1, sums1[[i]]$coefficients[2,2])
se2 = c(se2, sums2[[i]]$coefficients[2,2])
se3 = c(se3, sums3[[i]]$coefficients[2,2])
se4 = c(se4, sums4[[i]]$coefficients[2,2])
}
coef_df = data.frame("IV" = c(rep("Belongingness", length(dvs)),
rep("SelfEsteem", length(dvs)),
rep("Control", length(dvs)),
rep("MeaningfulExistence", length(dvs))),
"DV" = c(dvs,dvs,dvs,dvs),
"beta" = c(betas1, betas2,betas3,betas4),
"se" = c(se1, se2,se3,se4))
coef_df = coef_df
return(coef_df)
}dvs = c("subacc_roi_mask","rInsula_roi_mask",
"lInsula_roi_mask", "functional_combined_roi_mask")
coef_df = construct_coef_df_nts(dvs, df) %>%
mutate(DV = gsub("_roi_mask", "", DV))
plot_coef_iri(coef_dv)